# delimit ;
set more off ;
capture log close ;
clear all ;
version 16 ;

/* CODE FOR EXPERT ELICITATION ANALYSIS */
/* AUTHOR: ROGER H. VON HAEFEN  */
/* LAST EDITED: JANUARY 1, 2023 */

log using $ee_analysis_log, replace text ;
log off ;

/* EXPERT ELICITATION ANALYSIS */
/* EXPERT ELICITATION ANALYSIS */
/* EXPERT ELICITATION ANALYSIS */

/* ECOSYSTEM CONDITIONS & MURKY WATER DAYS */

use $ee_data, replace ;
drop if inlist(id,6,7,8) ;

gen ct = 3 ;
expand ct ;
drop ct ;
sort id question ;
by id question : gen cat = _n ;

gen double test = ec_good + ec_fair + ec_poor ;
assert test == 100 ;
replace ec_good = 100*ec_good/test ;
replace ec_fair = 100*ec_fair/test ;
replace ec_poor = 100*ec_poor/test ;
drop test ;

gen double test   = mw_low + mw_med + mw_high ;
sum test ;
replace mw_low    = 100*mw_low/test ;
replace mw_med    = 100*mw_med/test ;
replace mw_high   = 100*mw_high/test ;
drop test ;

gen ec = 0 ;
replace ec = round(ec_good) if cat == 1 ;
replace ec = round(ec_fair) if cat == 2 ;
replace ec = round(ec_poor) if cat == 3 ;

gen mw = 0 ;
replace mw = round(mw_low)    if cat == 1 ;
replace mw = round(mw_med)    if cat == 2 ;
replace mw = round(mw_high)   if cat == 3 ;

gen id1 = (id == 1) ; label var id1 "EXPERT #1" ;
gen id2 = (id == 2) ; label var id2 "EXPERT #2" ;
gen id3 = (id == 3) ; label var id3 "EXPERT #3" ;
gen id4 = (id == 4) ; label var id4 "EXPERT #4" ;
gen id5 = (id == 5) ; label var id5 "EXPERT #5" ;

order id question ec mw ;

/* ECOSYSTEM CONDITIONS */

foreach var in bi fc sc tn tp turb { ;
   sum `var' ;  replace `var' = (`var'-r(mean))/r(sd) ;
} ;   


preserve ; 
  drop if ec == 0 ;
  log on ;
  oprobit cat bi fc sc tn tp turb id1 id2 id3 id4 id5 [fw=ec], vce(cluster id) ;
  log off ;
  
  local bi_ec = _b[bi] ;
  local fc_ec = _b[fc] ;
  local sc_ec = _b[sc] ;
  local tn_ec = _b[tn] ;
  local tp_ec = _b[tp] ;
  local turb_ec = _b[turb] ;
  local id_ec = (_b[id1]+_b[id2]+_b[id3]+_b[id4]+_b[id5])/5 ;
  local cut1_ec = _b[/cut1] ;
  local cut2_ec = _b[/cut2] ;

  sum bi ; local bi_mean = r(mean) ; local bi_se = r(sd) ;
  sum fc ; local fc_mean = r(mean) ; local fc_se = r(sd) ;
  sum sc ; local sc_mean = r(mean) ; local sc_se = r(sd) ;
  sum tn ; local tn_mean = r(mean) ; local tn_se = r(sd) ;
  sum tp ; local tp_mean = r(mean) ; local tp_se = r(sd) ;
  sum turb ; local turb_mean = r(mean) ; local turb_se = r(sd) ;

  gen double index_ec0 = `bi_mean'*`bi_ec' + `fc_mean'*`fc_ec' +
                         `sc_mean'*`sc_ec' + `tn_mean'*`tn_ec' +
                         `tp_mean'*`tp_ec' + `id_ec' ;                     
  gen double cat1_ec0 = normal(`cut1_ec' - index_ec0) ;
  gen double cat2_ec0 = normal(`cut2_ec' - index_ec0) - normal(`cut1_ec' - index_ec0) ;
  gen double cat3_ec0 = 1 - normal(`cut2_ec' - index_ec0) ;
  
  local sub1 = "ec1" ;
  gen double index_`sub1' = index_ec0 + `bi_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;

  local sub1 = "ec2" ;
  gen double index_`sub1' = index_ec0 + `fc_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;

  local sub1 = "ec3" ;
  gen double index_`sub1' = index_ec0 + `sc_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;
  
  local sub1 = "ec4" ;
  gen double index_`sub1' = index_ec0 + `tn_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;

  local sub1 = "ec5" ;
  gen double index_`sub1' = index_ec0 + `tp_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;

  local sub1 = "ec6" ;
  gen double index_`sub1' = index_ec0 + `turb_ec'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_ec' - `ind') - cat1_ec0 ;
  gen double cat2_`sub1' = normal(`cut2_ec' - `ind') - normal(`cut1_ec' - `ind') - cat2_ec0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_ec' - `ind') - cat3_ec0 ;
  
  sum cat1* cat2* cat3* ;
  
restore ;

/* MURKY WATER DAYS */

preserve ;
  drop if mw == 0 ;
  log on ;
  oprobit cat bi fc sc tn tp turb id1 id2 id3 id4 id5 [fw=mw], vce(cluster id) ;
  log off ;

  local bi_mw = _b[bi] ;
  local fc_mw = _b[fc] ;
  local sc_mw = _b[sc] ;
  local tn_mw = _b[tn] ;
  local tp_mw = _b[tp] ;
  local turb_mw = _b[turb] ;
  local id_mw = (_b[id1]+_b[id2]+_b[id3]+_b[id4]+_b[id5])/5 ;
  local cut1_mw = _b[/cut1] ;
  local cut2_mw = _b[/cut2] ;

  sum bi ; local bi_mean = r(mean) ; local bi_se = r(sd) ;
  sum fc ; local fc_mean = r(mean) ; local fc_se = r(sd) ;
  sum sc ; local sc_mean = r(mean) ; local sc_se = r(sd) ;
  sum tn ; local tn_mean = r(mean) ; local tn_se = r(sd) ;
  sum tp ; local tp_mean = r(mean) ; local tp_se = r(sd) ;
  sum turb ; local turb_mean = r(mean) ; local turb_se = r(sd) ;

  gen double index_mw0 = `bi_mean'*`bi_mw' + `fc_mean'*`fc_mw' +
                         `sc_mean'*`sc_mw' + `tn_mean'*`tn_mw' +
                         `tp_mean'*`tp_mw' + `id_mw' ;                     
  gen double cat1_mw0 = normal(`cut1_mw' - index_mw0) ;
  gen double cat2_mw0 = normal(`cut2_mw' - index_mw0) - normal(`cut1_mw' - index_mw0) ;
  gen double cat3_mw0 = 1 - normal(`cut2_mw' - index_mw0) ;
  
  local sub1 = "mw1" ;
  gen double index_`sub1' = index_mw0 + `bi_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;

  local sub1 = "mw2" ;
  gen double index_`sub1' = index_mw0 + `fc_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;

  local sub1 = "mw3" ;
  gen double index_`sub1' = index_mw0 + `sc_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;
  
  local sub1 = "mw4" ;
  gen double index_`sub1' = index_mw0 + `tn_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;

  local sub1 = "mw5" ;
  gen double index_`sub1' = index_mw0 + `tp_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;

  local sub1 = "mw6" ;
  gen double index_`sub1' = index_mw0 + `turb_mw'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_mw' - `ind') - cat1_mw0 ;
  gen double cat2_`sub1' = normal(`cut2_mw' - `ind') - normal(`cut1_mw' - `ind') - cat2_mw0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_mw' - `ind') - cat3_mw0 ;
  
  sum cat1* cat2* cat3* ;

restore ;




/* HEALTH RISK */

use $ee_data, replace ;
sort question id ;
by question : assert bi == bi[1] ;
by question : assert fc == fc[1] ;
by question : assert sc == sc[1] ;
by question : assert tn == tn[1] ;
by question : assert tp == tp[1] ;
by question : assert turb == turb[1] ;

keep if inlist(id,6,7,8) ;

gen ct = 3 ;
expand ct ;
drop ct ;
sort id question ;
by id question : gen cat = _n ;

gen double test = hr_low + hr_medium + hr_high ;
assert test == 100 ;
drop test ;

gen hr = 0 ;
replace hr = round(hr_low)    if cat == 1 ;
replace hr = round(hr_medium) if cat == 2 ;
replace hr = round(hr_high)   if cat == 3 ;

order id question hr ;

gen id6 = (id == 6) ; label var id6 "EXPERT #6" ;
gen id7 = (id == 7) ; label var id7 "EXPERT #7" ;
gen id8 = (id == 8) ; label var id8 "EXPERT #8" ;

drop if hr == 0 ;

foreach var in bi fc sc tn tp turb { ;
   sum `var' ;  replace `var' = (`var'-r(mean))/r(sd) ;
} ;   

preserve ;
  log on ;
  oprobit cat bi  fc  sc  tn  tp  turb  id6 id7 id8 [fw=hr], vce(cluster id) ;  
  log off ;

  local bi_hr = _b[bi] ;
  local fc_hr = _b[fc] ;
  local sc_hr = _b[sc] ;
  local tn_hr = _b[tn] ;
  local tp_hr = _b[tp] ;
  local turb_hr = _b[turb] ;

  local id_hr = (_b[id6]+_b[id7]+_b[id8])/3 ;
  local cut1_hr = _b[/cut1] ;
  local cut2_hr = _b[/cut2] ;

  sum bi ; local bi_mean = r(mean) ; local bi_se = r(sd) ;
  sum fc ; local fc_mean = r(mean) ; local fc_se = r(sd) ;
  sum sc ; local sc_mean = r(mean) ; local sc_se = r(sd) ;
  sum tn ; local tn_mean = r(mean) ; local tn_se = r(sd) ;
  sum tp ; local tp_mean = r(mean) ; local tp_se = r(sd) ;
  sum turb ; local turb_mean = r(mean) ; local turb_se = r(sd) ;

  gen double index_hr0 = `bi_mean'*`bi_hr' + `fc_mean'*`fc_hr' +
                         `sc_mean'*`sc_hr' + `tn_mean'*`tn_hr' +
                         `tp_mean'*`tp_hr' + `id_hr' ;                     
  gen double cat1_hr0 = normal(`cut1_hr' - index_hr0) ;
  gen double cat2_hr0 = normal(`cut2_hr' - index_hr0) - normal(`cut1_hr' - index_hr0) ;
  gen double cat3_hr0 = 1 - normal(`cut2_hr' - index_hr0) ;
  
  local sub1 = "hr1" ;
  gen double index_`sub1' = index_hr0 + `bi_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;

  local sub1 = "hr2" ;
  gen double index_`sub1' = index_hr0 + `fc_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;

  local sub1 = "hr3" ;
  gen double index_`sub1' = index_hr0 + `sc_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;
  
  local sub1 = "hr4" ;
  gen double index_`sub1' = index_hr0 + `tn_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;

  local sub1 = "hr5" ;
  gen double index_`sub1' = index_hr0 + `tp_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;

  local sub1 = "hr6" ;
  gen double index_`sub1' = index_hr0 + `turb_hr'  ;                     
  local ind = "index_`sub1'" ;
  gen double cat1_`sub1' = normal(`cut1_hr' - `ind') - cat1_hr0 ;
  gen double cat2_`sub1' = normal(`cut2_hr' - `ind') - normal(`cut1_hr' - `ind') - cat2_hr0 ;
  gen double cat3_`sub1' = 1 - normal(`cut2_hr' - `ind') - cat3_hr0 ;
  
  sum cat1* cat2* cat3* ;

restore ;
